/**
* file name  : it/ciano/jenoware/SbjRelationship.java
* authors    : Luciano Xumerle
* created    : ven 02 giu 2006 17:22:50 CEST
*
* Copyright (c) 2004-2005 Luciano Xumerle. All rights reserved.
*
* This file is part of gephinfo.
*
* gephinfo is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program.  If not, see <http://www.gnu.org/licenses/>.
*
*/

package it.ciano.jenoware;

import java.util.HashMap;
import java.util.Arrays;
import it.ciano.util.*;

/**
 *
 * @author Luciano Xumerle
 * @version 0.0.1
 */
public class SbjRelationship extends HashMap
{

    static final long serialVersionUID = -1679388909033119496L;

    private boolean _DEBUG_;

    // store the MarkerFreq object
    CCounter mf;

    /**
     * Class constructor.
     *
     * @param alleleCounter the alleles counter.
     */
    public SbjRelationship ( CCounter alleleCounter )
    {
        this ( alleleCounter, false );
    }

    /**
     * Class constructor.
     *
     * @param alleleCounter the alleles counter.
     * @param useDebug set true to print debug information.
     */
    public SbjRelationship ( CCounter alleleCounter, boolean useDebug )
    {
        mf = alleleCounter;
        _DEBUG_ = useDebug;
    }

    /**
     * Class constructor.
     * The constructor uses the getUserAllelesCounter() method to get the CCounter.
     *
     * @param markerFreq the initialized MarkerFreq object.
     */
    public SbjRelationship ( MarkerFreq markerFreq )
    {
        this ( markerFreq.getUserAllelesCounter(), false );
    }

    /**
     * Class constructor.
     * The constructor uses the getUserAllelesCounter() method to get the CCounter.
     *
     * @param markerFreq the initialized MarkerFreq object.
     */
    public SbjRelationship ( MarkerFreq markerFreq, boolean useDebug )
    {
        this ( markerFreq.getUserAllelesCounter(), useDebug );
    }

    /**
     * Returns the probability of Genotype combinations for a Parent-Child relationships.
     *
     * @param a the PairAllele object.
     * @param b the PairAllele object.
     * @return the probability.
     */
    public double getParentChildProb ( PairAllele a, PairAllele b )
    {
        double[] res = getProbsList ( a, b );
        return res[ 0 ];
    }

    /**
     * Returns the probability of Genotype combinations for a Full Sibs relationships.
     *
     * @param a the PairAllele object.
     * @param b the PairAllele object.
     * @return the probability.
     */
    public double getFullSibsProb ( PairAllele a, PairAllele b )
    {
        double[] res = getProbsList ( a, b );
        return res[ 1 ];
    }

    /**
     * Returns the probability of Genotype combinations for a second degree relationships.
     *
     * @param a the PairAllele object.
     * @param b the PairAllele object.
     * @return the probability.
     */
    public double getSecondDegreeProb ( PairAllele a, PairAllele b )
    {
        double[] res = getProbsList ( a, b );
        return res[ 2 ];
    }

    /**
     * Returns the probability of Genotype combinations for a Non-relatives relationships.
     *
     * @param a the PairAllele object.
     * @param b the PairAllele object.
     * @return the probability.
     */
    public double getNonRelativesProb ( PairAllele a, PairAllele b )
    {
        double[] res = getProbsList ( a, b );
        return res[ 3 ];
    }

    final private double[] getProbsList ( PairAllele a, PairAllele b )
    {
        if ( a == null || b == null || !a.isValid() || !b.isValid() )
            return new double[] { -1, -1, -1, -1};

        PairAllele[] gg = new PairAllele[] { a, b };
        Arrays.sort ( gg );
        Allele[] aa = gg[0].getSortedAlleles();
        Allele[] bb = gg[1].getSortedAlleles();
        final String key = getGenotypesKey ( aa, bb );

        if ( !this.containsKey ( key ) )
        {
            // print debug info
            if ( _DEBUG_ )
                System.err.println ( "genotype key: " + key );
            // set the probs
            setGenotypesCat ( key, aa, bb );
        }
        return ( double[] ) this.get ( key );
    }


    final private void setGenotypesCat ( String key, Allele[] a, Allele[] b )
    {
        String a1;
        String a2;
        String b1;
        String b2;
        if (  b[0].equals ( b[1] ) && ! a[0].equals ( a[1] )  )
        {
            a1=b[0].toString();
            a2=b[1].toString();
            b1=a[0].toString();
            b2=a[1].toString();
        }
        else
        {
            a1=a[0].toString();
            a2=a[1].toString();
            b1=b[0].toString();
            b2=b[1].toString();
        }

        double g11 = mf.getFrequencies ( a1 );
        double g12 = mf.getFrequencies ( a2 );
        double g21 = mf.getFrequencies ( b1 );
        double g22 = mf.getFrequencies ( b2 );

        double res[];

        // the first is homozigote
        if ( a1.equals ( a2 ) )
        {
            // the second is homozigote
            if ( b1.equals ( b2 ) )
            {
                if ( a1.equals ( b1 ) )
                {
                    res = calculateProbs ( 1, g11, g12, g21, g22 );
                }
                else
                {
                    res = calculateProbs ( 3, g11, g12, g21, g22 );
                }
            }
            else
            {
                if ( a1.equals ( b1 ) )
                {
                    res = calculateProbs ( 2, g11, g12, g21, g22 );
                }
                else if ( a1.equals ( b2 ) )
                {
                    res = calculateProbs ( 2, g11, g12, g22, g21 );
                }
                else
                {
                    res = calculateProbs ( 5, g11, g12, g21, g22 );
                }
            }
        }
        else
        {
            if ( a1.equals ( b1 ) && a2.equals ( b2 ) )
            {
                res = calculateProbs ( 4, g11, g12, g21, g22 );
            }
            else if ( a1.equals ( b1 ) )
            {
                res = calculateProbs ( 6, g11, g12, g21, g22 );
            }
            else if ( a1.equals ( b2 ) )
            {
                res = calculateProbs ( 6, g11, g12, g22, g21 );
            }
            else if ( a2.equals ( b1 ) )
            {
                res = calculateProbs ( 6, g12, g11, g21, g22 );
            }
            else if ( a2.equals ( b2 ) )
            {
                res = calculateProbs ( 6, g12, g11, g22, g21 );
            }
            else
            {
                res = calculateProbs ( 7, g11, g12, g21, g22 );
            }
        }
        this.put ( key, res );
    }


    final private static String getGenotypesKey ( Allele[] aa, Allele[] bb )
    {
        StringBuffer res = new StringBuffer();
        res.append ( "-" );
        res.append ( aa[0].toString() );
        res.append ( "-" );
        res.append ( aa[1].toString() );
        res.append ( "-" );
        res.append ( bb[0].toString() );
        res.append ( "-" );
        res.append ( bb[1].toString() );
        res.append ( "-" );
        return res.toString();
    }

    /**
     * Returns an array containing 4 probs:
     * { parent-child, full_sibs, second_degree, non_relatives  }
     *
     * @param type the type of relation (see reference).
     * @param a allele freq.
     * @param b allele freq.
     * @param c allele freq.
     * @param d allele freq.
     * @return the array of probabilities.
     */
    final private double[] calculateProbs ( int type, double a, double b, double c, double d )
    {
        double res[] = new double[] { -1, -1, -1, -1};

        // print debug informations
        if ( _DEBUG_ )
            System.err.println ( " type: " + type + "\nGeno1 Freq: "
                    + a + " - " + b + "\nGeno2 Freq: " + c + " - " + d );

        switch ( type )
        {
            // AA,AA
        case 1:
            res[ 0 ] = a * b * c;
            res[ 1 ] = a * b * ( 1 + a ) * ( 1 + a ) / 4;
            res[ 2 ] = res[ 0 ] * ( 1 + a ) / 2;
            res[ 3 ] = res[ 0 ] * d;
            break;
            // AA,AB
        case 2:
            res[ 0 ] = 2 * a * b * d;
            res[ 1 ] = a * b * d * ( 1 + a );
            res[ 2 ] = a * b * d * ( 1 + 2 * a );
            res[ 3 ] = 4 * a * b * c * d;
            break;
            // AA,BB
        case 3:
            res[ 0 ] = 0;
            res[ 1 ] = a * b * c * d / 2;
            res[ 2 ] = a * b * c * d;
            res[ 3 ] = 2 * res[ 2 ];
            break;
            // AB,AB
        case 4:
            res[ 0 ] = a * b * ( a + b );
            res[ 1 ] = a * b * ( 2 * a * b + a + b + 1 ) / 2;
            res[ 2 ] = a * b * ( 4 * a * b + a + b ) / 2;
            res[ 3 ] = 4 * a * b * c * d;
            break;
            // AA,BC
        case 5:
            res[ 0 ] = 0;
            res[ 1 ] = a * b * c * d;
            res[ 2 ] = 2 * res[ 1 ];
            res[ 3 ] = 4 * res[ 1 ];
            break;
            // AB,AC
        case 6:
            res[ 0 ] = 2 * a * b * d;
            res[ 1 ] = a * b * d * ( 2 * a + 1 );
            res[ 2 ] = a * b * d * ( 4 * a + 1 );
            res[ 3 ] = 8 * a * b * c * d;
            break;
            // AB,CD
        case 7:
            res[ 0 ] = 0;
            res[ 1 ] = 2 * a * b * c * d;
            res[ 2 ] = 4 * a * b * c * d;
            res[ 3 ] = 8 * a * b * c * d;
            break;
        default:
            res = new double[]
            {
                -1, -1, -1, -1
            };
            break;
        }
        // debug results
        if ( _DEBUG_ )
            System.err.println (  "\np-c: "  + res[ 0 ] + "\nf-s: " + res[ 1 ]
                    + "\n2nd: " + res[ 2 ] + "\nunr: " + res[ 3 ] );
        return res;
    }




    /**
     * Overriding toString method.
     *
     * @return the String.
     */
    public String toString()
    {
        return "";
    }

} // end class

